Chapter 5.4

library(tidyverse) # metapackage of all tidyverse packages

## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.0 ──

## ✓ ggplot2 3.3.3     ✓ purrr   0.3.4
## ✓ tibble  3.0.4     ✓ dplyr   1.0.2
## ✓ tidyr   1.1.2     ✓ stringr 1.4.0
## ✓ readr   1.4.0     ✓ forcats 0.5.0

## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## x dplyr::filter() masks stats::filter()
## x dplyr::lag()    masks stats::lag()

library(keras)
use_condaenv("r-reticulate")

5.4.1

Load our model

model <- load_model_hdf5("cats_and_dogs_small_2.h5")
summary(model)

## Model: "sequential_2"
## ________________________________________________________________________________
## Layer (type)                        Output Shape                    Param #     
## ================================================================================
## conv2d_10 (Conv2D)                  (None, 148, 148, 32)            896         
## ________________________________________________________________________________
## max_pooling2d_9 (MaxPooling2D)      (None, 74, 74, 32)              0           
## ________________________________________________________________________________
## conv2d_9 (Conv2D)                   (None, 72, 72, 64)              18496       
## ________________________________________________________________________________
## max_pooling2d_8 (MaxPooling2D)      (None, 36, 36, 64)              0           
## ________________________________________________________________________________
## conv2d_8 (Conv2D)                   (None, 34, 34, 128)             73856       
## ________________________________________________________________________________
## max_pooling2d_7 (MaxPooling2D)      (None, 17, 17, 128)             0           
## ________________________________________________________________________________
## conv2d_7 (Conv2D)                   (None, 15, 15, 128)             147584      
## ________________________________________________________________________________
## max_pooling2d_6 (MaxPooling2D)      (None, 7, 7, 128)               0           
## ________________________________________________________________________________
## flatten_2 (Flatten)                 (None, 6272)                    0           
## ________________________________________________________________________________
## dropout (Dropout)                   (None, 6272)                    0           
## ________________________________________________________________________________
## dense_5 (Dense)                     (None, 512)                     3211776     
## ________________________________________________________________________________
## dense_4 (Dense)                     (None, 1)                       513         
## ================================================================================
## Total params: 3,453,121
## Trainable params: 3,453,121
## Non-trainable params: 0
## ________________________________________________________________________________

Get input image:

img_path <- "cats_and_dogs_small/test/cats/cat.1700.jpg"
img <- image_load(img_path, target_size = c(150, 150))                 
img_tensor <- image_to_array(img)
img_tensor <- array_reshape(img_tensor, c(1, 150, 150, 3))
img_tensor <- img_tensor / 255                                         
dim(img_tensor)                                                        

## [1]   1 150 150   3

take a look at the image

plot(as.raster(img_tensor[1,,,]))

now create the model. Using keras_model instead of keras_sequential_model allows us to access multiple output layers

layer_outputs <- lapply(model$layers[1:8], function(layer) layer$output)      
activation_model <- keras_model(inputs = model$input, outputs = layer_outputs)

str(layer_outputs) # at this point I think it is just a list of (empty) tensors

## List of 8
##  $ :Tensor("conv2d_10/Identity:0", shape=(None, 148, 148, 32), dtype=float32)
##  $ :Tensor("max_pooling2d_9/Identity:0", shape=(None, 74, 74, 32), dtype=float32)
##  $ :Tensor("conv2d_9/Identity:0", shape=(None, 72, 72, 64), dtype=float32)
##  $ :Tensor("max_pooling2d_8/Identity:0", shape=(None, 36, 36, 64), dtype=float32)
##  $ :Tensor("conv2d_8/Identity:0", shape=(None, 34, 34, 128), dtype=float32)
##  $ :Tensor("max_pooling2d_7/Identity:0", shape=(None, 17, 17, 128), dtype=float32)
##  $ :Tensor("conv2d_7/Identity:0", shape=(None, 15, 15, 128), dtype=float32)
##  $ :Tensor("max_pooling2d_6/Identity:0", shape=(None, 7, 7, 128), dtype=float32)

activations <- activation_model %>% predict(img_tensor)         

str(activations)

## List of 8
##  $ : num [1, 1:148, 1:148, 1:32] 0.00837 0.01414 0.00918 0.01291 0.00825 ...
##  $ : num [1, 1:74, 1:74, 1:32] 0.0141 0.0133 0.0178 0.021 0.0205 ...
##  $ : num [1, 1:72, 1:72, 1:64] 0 0 0 0 0 0 0 0 0 0 ...
##  $ : num [1, 1:36, 1:36, 1:64] 0 0 0 0 0 0 0 0 0 0 ...
##  $ : num [1, 1:34, 1:34, 1:128] 0.00386 0 0 0 0 ...
##  $ : num [1, 1:17, 1:17, 1:128] 0.00386 0 0 0.0131 0.00835 ...
##  $ : num [1, 1:15, 1:15, 1:128] 0.01514 0 0.01075 0.00422 0.01076 ...
##  $ : num [1, 1:7, 1:7, 1:128] 0.02181 0.01075 0.01471 0.00678 0.02665 ...

define plotting function

plot_channel <- function(channel) {
  rotate <- function(x) t(apply(x, 2, rev))
  image(rotate(channel), axes = FALSE, asp = 1,
        col = topo.colors(12))
}

first_layer_activation <- activations[[1]]
dim(first_layer_activation)

## [1]   1 148 148  32

plot_channel(first_layer_activation[1,,,2])

plot_channel(first_layer_activation[1,,,7])

plot them all

image_size <- 58
images_per_row <- 16

for (i in 1:8) {

  layer_activation <- activations[[i]]
  layer_name <- model$layers[[i]]$name

  n_features <- dim(layer_activation)[[4]]
  n_cols <- n_features %/% images_per_row

  #png(paste0("cat_activations_", i, "_", layer_name, ".png"),
   #   width = image_size * images_per_row,
  #    height = image_size * n_cols)
  op <- par(mfrow = c(n_cols, images_per_row), mai = rep_len(0.02, 4))

  for (col in 0:(n_cols-1)) {
    for (row in 0:(images_per_row-1)) {
      channel_image <- layer_activation[1,,,(col*images_per_row) + row + 1]
      plot_channel(channel_image)
    }
  }

  par(op)
  #dev.off()
}

5.4.2

visualizing the filters

set up the loss function

library(keras)
library(tensorflow)
tf$compat$v1$disable_eager_execution()
model <- application_vgg16(
  weights = "imagenet",
  include_top = FALSE
)
layer_name <- "block3_conv1"
filter_index <- 1
layer_output <- get_layer(model, layer_name)$output
loss <- k_mean(layer_output[,,,filter_index]) # average output as a tensor

get the gradient associated with the above loss and normalize (RMS)

grads <- k_gradients(loss, model$input)[[1]]  
grads <- grads / (k_sqrt(k_mean(k_square(grads))) + 1e-5) # as a tensor

now we need to be able to calculate loss and gradient for a given input. We use iterate for this:

iterate <- k_function(list(model$input), list(loss, grads))
c(loss_value, grads_value) %<-%
    iterate(list(array(0, dim = c(1, 150, 150, 3))))

put it together into a loop

input_img_data <-                                                     
  array(runif(150 * 150 * 3), dim = c(1, 150, 150, 3)) * 20 + 128    # input random image, near grey 
step <- 1
for (i in 1:40) {                                                     
  c(loss_value, grads_value) %<-% iterate(list(input_img_data)) # calculate gradient and loss
  cat("loss: ", loss_value, "\n")
  cat("grads_value: ", grads_value[1,1:5,1,1], "\n")
  input_img_data <- input_img_data + (grads_value * step) # update image     
}

## loss:  71.47847 
## grads_value:  -0.0359814 0.004708725 0.07507561 0.03921726 0.03366342 
## loss:  154.1098 
## grads_value:  -0.03196237 0.002609856 0.06951085 0.04309835 0.03370903 
## loss:  257.8853 
## grads_value:  -0.05768544 0.03202096 0.1526623 0.1001099 0.04960227 
## loss:  380.5151 
## grads_value:  -0.05373624 0.04652522 0.1467815 0.07892174 0.05268444 
## loss:  514.2333 
## grads_value:  -0.06410471 0.03976501 0.1656701 0.09998608 0.06395946 
## loss:  651.1523 
## grads_value:  -0.05772654 0.04214632 0.1538958 0.1029205 0.09177531 
## loss:  788.8897 
## grads_value:  -0.06127032 0.04556536 0.1535099 0.08913926 0.09020352 
## loss:  925.9111 
## grads_value:  -0.04143165 0.05982845 0.1656604 0.1185635 0.07577652 
## loss:  1061.323 
## grads_value:  -0.04561591 0.05401523 0.1648132 0.1299802 0.08686046 
## loss:  1195.649 
## grads_value:  -0.03925656 0.06155571 0.177838 0.1550038 0.09009928 
## loss:  1329.272 
## grads_value:  -0.03801572 0.04789129 0.1378204 0.1362879 0.06986129 
## loss:  1462.021 
## grads_value:  -0.04677001 0.0569244 0.1596694 0.1536451 0.07338627 
## loss:  1593.7 
## grads_value:  -0.04448375 0.04576671 0.1492003 0.1327797 0.05765094 
## loss:  1724.606 
## grads_value:  -0.05825286 0.01009258 0.1185173 0.1371254 0.06848843 
## loss:  1854.169 
## grads_value:  -0.05572954 0.0101013 0.1151108 0.1209305 0.04711407 
## loss:  1982.357 
## grads_value:  -0.03328842 0.03709821 0.1556707 0.1158074 0.04246084 
## loss:  2108.965 
## grads_value:  -0.01265447 0.05834574 0.1631498 0.08730896 0.02336429 
## loss:  2234.811 
## grads_value:  -0.02815587 0.04371713 0.1665642 0.1030919 0.05894537 
## loss:  2359.484 
## grads_value:  -0.01330472 0.08740627 0.22189 0.1485916 0.1039804 
## loss:  2483.485 
## grads_value:  -0.01879269 0.0601967 0.2478185 0.1913538 0.1923828 
## loss:  2606.848 
## grads_value:  0.02221461 0.1329465 0.3152845 0.1275125 0.07358015 
## loss:  2729.285 
## grads_value:  0.04912999 0.1051024 0.2981746 -0.02365829 -0.04870183 
## loss:  2851.139 
## grads_value:  0.074811 0.1551519 0.3499207 0.09673315 0.06032985 
## loss:  2972.572 
## grads_value:  0.08059686 0.1686005 0.3473207 -0.1256245 -0.1747828 
## loss:  3093.191 
## grads_value:  0.05231205 0.2076476 0.4489791 0.05413037 -0.04625894 
## loss:  3213.107 
## grads_value:  0.07596014 0.2366392 0.4494065 0.03994574 -0.05891027 
## loss:  3332.267 
## grads_value:  0.08621653 0.2447456 0.4377641 0.02927187 -0.07044771 
## loss:  3450.695 
## grads_value:  0.103541 0.2337443 0.4260218 0.02727507 -0.1018605 
## loss:  3568.775 
## grads_value:  0.04630553 0.1526673 0.3619846 -0.1306613 -0.2013386 
## loss:  3686.101 
## grads_value:  0.04027896 0.1564929 0.3832963 0.02326949 -0.08941849 
## loss:  3802.767 
## grads_value:  0.06704955 0.1272947 0.3403085 -0.2216826 -0.2566642 
## loss:  3918.843 
## grads_value:  0.06692022 0.1595799 0.3839905 -0.107966 -0.1620809 
## loss:  4034.278 
## grads_value:  0.06720355 0.1876347 0.4107423 -0.1082291 -0.1894303 
## loss:  4149.093 
## grads_value:  0.05971123 0.1815435 0.376109 -0.03583135 -0.153118 
## loss:  4263.442 
## grads_value:  0.06043259 0.1743895 0.3747339 -0.04567286 -0.11435 
## loss:  4377.311 
## grads_value:  0.03362913 0.1413212 0.3807378 -0.04882617 -0.1557696 
## loss:  4490.654 
## grads_value:  0.02378552 0.1452153 0.3896282 -0.005004994 -0.2355956 
## loss:  4603.606 
## grads_value:  0.02875094 0.1564975 0.4024472 0.009540815 -0.2339753 
## loss:  4716.061 
## grads_value:  0.03596384 0.1222304 0.3634512 0.05770922 -0.2513396 
## loss:  4828.246 
## grads_value:  0.03012265 0.1111424 0.3485721 0.0907127 -0.2067998

gradient ascent because we are increasing the loss?

post process the tensor so that we can dispaly it as an image:

deprocess_image <- function(x) {
  dms <- dim(x)
  x <- x - mean(x)                        
  x <- x / (sd(x) + 1e-5)                 
  x <- x * 0.1                            
  x <- x + 0.5                            
  x <- pmax(0, pmin(x, 1))                
  array(x, dim = dms)                     
}

put it all together in a function

generate_pattern <- function(layer_name, filter_index, size = 150) {
  layer_output <- model$get_layer(layer_name)$output                      
  loss <- k_mean(layer_output[,,,filter_index])                           
  grads <- k_gradients(loss, model$input)[[1]]                            
  grads <- grads / (k_sqrt(k_mean(k_square(grads))) + 1e-5)               
  iterate <- k_function(list(model$input), list(loss, grads))             
  input_img_data <-                                                       
    array(runif(size * size * 3), dim = c(1, size, size, 3)) * 20 + 128   
  step <- 1                                                               
  for (i in 1:40) {                                                       
    c(loss_value, grads_value) %<-% iterate(list(input_img_data))         
    input_img_data <- input_img_data + (grads_value * step)               
  }                                                                       
  img <- input_img_data[1,,,]
  deprocess_image(img)
}

library(grid)
grid.raster(generate_pattern("block3_conv1", 1))

library(grid)
library(gridExtra)
dir.create("vgg_filters")
for (layer_name in c("block1_conv1", "block2_conv1",
                     "block3_conv1", "block4_conv1")) {
  size <- 140

  png(paste0("vgg_filters/", layer_name, ".png"),
      width = 8 * size, height = 8 * size)

  grobs <- list()
  for (i in 0:7) {
    for (j in 0:7) {
      pattern <- generate_pattern(layer_name, i + (j*8) + 1, size = size)
      grob <- rasterGrob(pattern,
                         width = unit(0.9, "npc"),
                         height = unit(0.9, "npc"))
      grobs[[length(grobs)+1]] <- grob
    }
  }

  grid.arrange(grobs = grobs, ncol = 8)
  dev.off()
}

5.4.3 heat maps

model <- application_vgg16(weights = "imagenet")            

img_path <- "creative_commons_elephant.jpg"              
img <- image_load(img_path, target_size = c(224, 224)) %>%           
  image_to_array() %>%                                               
  array_reshape(dim = c(1, 224, 224, 3)) %>%                         
  imagenet_preprocess_input()                                        

 preds <- model %>% predict(img)

 imagenet_decode_predictions(preds, top = 3)[[1]]

which.max(preds[1,])

## [1] 387

summary(model)

## Model: "vgg16"
## ________________________________________________________________________________
## Layer (type)                        Output Shape                    Param #     
## ================================================================================
## input_2 (InputLayer)                [(None, 224, 224, 3)]           0           
## ________________________________________________________________________________
## block1_conv1 (Conv2D)               (None, 224, 224, 64)            1792        
## ________________________________________________________________________________
## block1_conv2 (Conv2D)               (None, 224, 224, 64)            36928       
## ________________________________________________________________________________
## block1_pool (MaxPooling2D)          (None, 112, 112, 64)            0           
## ________________________________________________________________________________
## block2_conv1 (Conv2D)               (None, 112, 112, 128)           73856       
## ________________________________________________________________________________
## block2_conv2 (Conv2D)               (None, 112, 112, 128)           147584      
## ________________________________________________________________________________
## block2_pool (MaxPooling2D)          (None, 56, 56, 128)             0           
## ________________________________________________________________________________
## block3_conv1 (Conv2D)               (None, 56, 56, 256)             295168      
## ________________________________________________________________________________
## block3_conv2 (Conv2D)               (None, 56, 56, 256)             590080      
## ________________________________________________________________________________
## block3_conv3 (Conv2D)               (None, 56, 56, 256)             590080      
## ________________________________________________________________________________
## block3_pool (MaxPooling2D)          (None, 28, 28, 256)             0           
## ________________________________________________________________________________
## block4_conv1 (Conv2D)               (None, 28, 28, 512)             1180160     
## ________________________________________________________________________________
## block4_conv2 (Conv2D)               (None, 28, 28, 512)             2359808     
## ________________________________________________________________________________
## block4_conv3 (Conv2D)               (None, 28, 28, 512)             2359808     
## ________________________________________________________________________________
## block4_pool (MaxPooling2D)          (None, 14, 14, 512)             0           
## ________________________________________________________________________________
## block5_conv1 (Conv2D)               (None, 14, 14, 512)             2359808     
## ________________________________________________________________________________
## block5_conv2 (Conv2D)               (None, 14, 14, 512)             2359808     
## ________________________________________________________________________________
## block5_conv3 (Conv2D)               (None, 14, 14, 512)             2359808     
## ________________________________________________________________________________
## block5_pool (MaxPooling2D)          (None, 7, 7, 512)               0           
## ________________________________________________________________________________
## flatten (Flatten)                   (None, 25088)                   0           
## ________________________________________________________________________________
## fc1 (Dense)                         (None, 4096)                    102764544   
## ________________________________________________________________________________
## fc2 (Dense)                         (None, 4096)                    16781312    
## ________________________________________________________________________________
## predictions (Dense)                 (None, 1000)                    4097000     
## ================================================================================
## Total params: 138,357,544
## Trainable params: 138,357,544
## Non-trainable params: 0
## ________________________________________________________________________________

african_elephant_output <- model$output[, 387]                             
last_conv_layer <- model %>% get_layer("block5_conv3")                     
grads <- k_gradients(african_elephant_output, last_conv_layer$output)[[1]] 
pooled_grads <- k_mean(grads, axis = c(1, 2, 3))                           
iterate <- k_function(list(model$input),                                   
                      list(pooled_grads, last_conv_layer$output[1,,,]))
c(pooled_grads_value, conv_layer_output_value) %<-% iterate(list(img))     
for (i in 1:512) {       # 512 channels                                             
  conv_layer_output_value[,,i] <-
    conv_layer_output_value[,,i] * pooled_grads_value[[i]]
}
heatmap <- apply(conv_layer_output_value, c(1,2), mean)                    

heatmap <- pmax(heatmap, 0)
heatmap <- heatmap / max(heatmap)                                          
write_heatmap <- function(heatmap, filename, width = 224, height = 224,    
                          bg = "white", col = terrain.colors(12)) {
  png(filename, width = width, height = height, bg = bg)
  op = par(mar = c(0,0,0,0))
  on.exit({par(op); dev.off()}, add = TRUE)
  rotate <- function(x) t(apply(x, 2, rev))
  image(rotate(heatmap), axes = FALSE, asp = 1, col = col)
}
write_heatmap(heatmap, "elephant_heatmap.png")                             

library(magick)

## Linking to ImageMagick 6.9.11.32
## Enabled features: cairo, fontconfig, freetype, lcms, pango, rsvg, webp
## Disabled features: fftw, ghostscript, x11

library(viridis)

## Loading required package: viridisLite

image <- image_read(img_path)                                      
info <- image_info(image)
geometry <- sprintf("%dx%d!", info$width, info$height)
pal <- col2rgb(viridis(20), alpha = TRUE)                          
alpha <- floor(seq(0, 255, length = ncol(pal)))
pal_col <- rgb(t(pal), alpha = alpha, maxColorValue = 255)
write_heatmap(heatmap, "elephant_overlay.png",
              width = 14, height = 14, bg = NA, col = pal_col)
image_read("elephant_overlay.png") %>%                             
  image_resize(geometry, filter = "quadratic") %>%
  image_composite(image, operator = "blend", compose_args = "20") %>%
  plot()